// Load the data.
use "$pathD/all_tract.dta", clear

// Generate dummies for income quartile.
tab income_quart, gen(income_quart)

// Generate and label the outoor PM income interactions.
gen med_income_log_outdoor = above_med_income * log_PM_2p5_outdoor
label var med_income_log_outdoor ///
	"Log outdoor PM\textsubscript{2.5} $\times$ above-median income"
gen med_income_outdoor = above_med_income * PM_2p5_outdoor
label var med_income_outdoor ///
	"Outdoor PM\textsubscript{2.5} $\times$ above-median income"
forvalues i = 2/4 {
	gen income_quart`i'_log_outdoor = income_quart`i' * log_PM_2p5_outdoor
	label var income_quart`i'_log_outdoor ///
		"Log outdoor PM\textsubscript{2.5} $\times$ `i'.\ income quartile"
	gen income_quart`i'_outdoor = income_quart`i' * PM_2p5_outdoor
	label var income_quart`i'_outdoor ///
		"Log outdoor PM\textsubscript{2.5} $\times$ `i'.\ income quartile"
	}
gen log_income_log_outoor = log_med_hh_income * log_PM_2p5_outdoor
label var log_income_log_outoor ///
	"Log outdoor PM\textsubscript{2.5} $\times$ log income"
gen log_income_outoor = log_med_hh_income * PM_2p5_outdoor
label var log_income_outoor "Outdoor PM\textsubscript{2.5} $\times$ log income"
	
// Generate dummies for wildfire smoke and interactions with income.
forvalues i = 1/2 {
	gen smokepm_pred`i' = med_smokepm_pred == `i'
	gen med_income_smoke`i' = smokepm_pred`i' * above_med_income
	gen log_income_smoke`i' = smokepm_pred`i' * log_med_hh_income
	}
label var smokepm_pred1 "Below-median wildfire PM\textsubscript{2.5}"
label var smokepm_pred2 "Above-median wildfire PM\textsubscript{2.5}"
label var med_income_smoke1 ///
	"Below-median wildfire PM\textsubscript{2.5} $\times$ above-median income"
label var med_income_smoke2 ///
	"Above-median wildfire PM\textsubscript{2.5} $\times$ above-median income"
label var log_income_smoke1 ///
	"Below-median wildfire PM\textsubscript{2.5} $\times$ log income"
label var log_income_smoke2 ///
	"Above-median wildfire PM\textsubscript{2.5} $\times$ log income"
forvalues i = 2/4 {
	gen income_quart`i'_smoke1 = income_quart`i' * smokepm_pred1
	label var income_quart`i'_smoke1 ///
		"Below-median wildfire PM\textsubscript{2.5} $\times$ `i'.\ income quartile"
	}
forvalues i = 2/4 {
	gen income_quart`i'_smoke2 = income_quart`i' * smokepm_pred2
	label var income_quart`i'_smoke2 ///
		"Above-median wildfire PM\textsubscript{2.5} $\times$ `i'.\ income quartile"
	}

// Generate the temporal fixed effects.
gen month = month(date)
gen dow = dow(date)

// Set the globals.
global weather_covs = ///
	"AIR_temp_* AIR_dew_* PRECIP_quant* WIND_speed_*"
global time = "i.month i.dow"
global opt = ", fe cluster(tract_id)"

// Set the panel variable.
destring tract_id, replace
xtset tract_id

// Main regressions
eststo: xtreg log_PM_2p5_indoor log_PM_2p5_outdoor med_income_log_outdoor ///
	$weather_covs $time $opt
eststo: xtreg log_PM_2p5_indoor log_PM_2p5_outdoor log_income_log_outoor ///
	$weather_covs $time $opt
eststo: xtreg log_PM_2p5_indoor smokepm_pred1 smokepm_pred2 ///
	med_income_smoke1 med_income_smoke2 ///
	$weather_covs $time $opt
eststo: xtreg log_PM_2p5_indoor smokepm_pred1 smokepm_pred2 ///
	log_income_smoke1 log_income_smoke2 ///
	$weather_covs $time $opt

// Table outdoor regressions
estout est1 est2 est3 est4 using ///
	"$pathR/Tables/main.tex", replace ///
	keep(log_PM_2p5_outdoor med_income_log_outdoor ///
	log_income_log_outoor smokepm_pred1 smokepm_pred2 ///
	med_income_smoke1 med_income_smoke2 ///
	log_income_smoke1 log_income_smoke2) ///
	order() ///
	label starl(* 0.10 ** 0.05 *** 0.01) mlabels(, prefix(%) begin(%)) ///
	varlabels(, end(" \addlinespace[0pt]")) ///
	cells(b(star fmt(3)) se(fmt(3) par)) style(tex) collabels(none) ///
	eqlabels(none)
estout est1 est2 est3 est4 using ////
	"$pathR/Tables/main_N.tex", replace ///
	keep(log_PM_2p5_outdoor) label omitted starl(* 0.10 ** 0.05 *** 0.01) ///
	mlabels(, title prefix(%) begin(%)) varlabels(, prefix(%) end()) ///
	cells() style(tex) collabels(none) stats(N, layout(@ "{@}") ///
	fmt(%12.0gc) labels("Observations")) eqlabels(none)
eststo clear
	
// Quartile regressions
eststo: xtreg log_PM_2p5_indoor log_PM_2p5_outdoor income_quart2_log_outdoor ///
	income_quart3_log_outdoor income_quart4_log_outdoor ///
	$weather_covs $time $opt
eststo: xtreg log_PM_2p5_indoor smokepm_pred1 smokepm_pred2 ///
	income_quart2_smoke1 income_quart3_smoke1 income_quart4_smoke1 ///
	income_quart2_smoke2 income_quart3_smoke2 income_quart4_smoke2 ///
	$weather_covs $time $opt

// Table quartile regressions
estout est1 est2 using ///
	"$pathR/Tables/quartile.tex", replace ///
	keep(log_PM_2p5_outdoor income_quart2_log_outdoor ///
	income_quart3_log_outdoor income_quart4_log_outdoor ///
	smokepm_pred1 smokepm_pred2 ///
	income_quart2_smoke1 income_quart3_smoke1 income_quart4_smoke1 ///
	income_quart2_smoke2 income_quart3_smoke2 income_quart4_smoke2) ///
	order() ///
	label starl(* 0.10 ** 0.05 *** 0.01) mlabels(, prefix(%) begin(%)) ///
	varlabels(, end(" \addlinespace[0pt]")) ///
	cells(b(star fmt(3)) se(fmt(3) par)) style(tex) collabels(none) ///
	eqlabels(none)
estout est1 est2 using ////
	"$pathR/Tables/quartile_N.tex", replace ///
	keep(log_PM_2p5_outdoor) label omitted starl(* 0.10 ** 0.05 *** 0.01) ///
	mlabels(, title prefix(%) begin(%)) varlabels(, prefix(%) end()) ///
	cells() style(tex) collabels(none) stats(N, layout(@ "{@}") ///
	fmt(%12.0gc) labels("Observations")) eqlabels(none)
eststo clear

// Linear regressions
eststo: xtreg PM_2p5_indoor PM_2p5_outdoor med_income_outdoor ///
	$weather_covs $time $opt
eststo: xtreg PM_2p5_indoor PM_2p5_outdoor log_income_outoor ///
	$weather_covs $time $opt
eststo: xtreg PM_2p5_indoor smokepm_pred1 smokepm_pred2 ///
	med_income_smoke1 med_income_smoke2 ///
	$weather_covs $time $opt
eststo: xtreg PM_2p5_indoor smokepm_pred1 smokepm_pred2 ///
	log_income_smoke1 log_income_smoke2 ///
	$weather_covs $time $opt
	
// Table linear regressions
estout est1 est2 est3 est4 using ///
	"$pathR/Tables/linear.tex", replace ///
	keep(PM_2p5_outdoor med_income_outdoor log_income_outoor ///
	smokepm_pred1 smokepm_pred2 ///
	med_income_smoke1 med_income_smoke2 ///
	log_income_smoke1 log_income_smoke2) ///
	order() ///
	label starl(* 0.10 ** 0.05 *** 0.01) mlabels(, prefix(%) begin(%)) ///
	varlabels(, end(" \addlinespace[0pt]")) ///
	cells(b(star fmt(3)) se(fmt(3) par)) style(tex) collabels(none) ///
	eqlabels(none)
estout est1 est2 est3 est4 using ////
	"$pathR/Tables/linear_N.tex", replace ///
	keep(PM_2p5_outdoor) label omitted starl(* 0.10 ** 0.05 *** 0.01) ///
	mlabels(, title prefix(%) begin(%)) varlabels(, prefix(%) end()) ///
	cells() style(tex) collabels(none) stats(N, layout(@ "{@}") ///
	fmt(%12.0gc) labels("Observations")) eqlabels(none)
eststo clear

clear
